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ABSTRACT 

We use the extended Press-Schechter formalism to investigate the rate at which cold 
dark matter haloes accrete mass. We discuss the shortcomings of previous methods 
that have been used to compute the mass accretion histories of dark matter haloes, 
and present an improved method based on the ./V-branch merger tree algorithm of 
Somerville & Kolatt. We show that this method no longer suffers from inconsistencies 
in halo formation times, and compare its predictions with high resolution iV-body sim- 
ulations. Although the overall agreement is reasonable, there are slight inconsistencies 
which are most easily interpreted as a reflection of ellipsoidal collapse (as opposed to 
spherical collapse assumed in the Press-Schechter formalism). We show that the aver- 
age mass accretion histories follow a simple, universal profile, and we present a simple 
recipe for computing the two scale-parameters which is applicable to a wide range of 
halo masses and cosmologies. Together with the universal profiles for the density and 
angular momentum distributions of CDM haloes, these universal mass accretion histo- 
ries provide a simple but accurate framework for modeling the structure and formation 
of dark matter haloes. In particular, they can be used as a backbone for modeling var- 
ious aspects of galaxy formation where one is not interested in the detailed effects of 
merging. As an example we use the universal mass accretion history to compute the 
rate at which dark matter haloes accrete mass, which we compare to the cosmic star 
formation history of the Universe. 

Key words: cosmology: theory — galaxies: formation — galaxies: halos — stars: 
formation — dark matter. 



1 INTRODUCTION 

In the standard Cold Dark Matter (CDM) family of cos- 
mological models dark matter haloes form hierarchically 
through the accretion and merging of smaller structures that 
condense out of a Gaussian initial density field. The rate at 
which these dark matter haloes grow in mass sets, amongst 
others, the rate at which baryons can cool to form luminous 
objects. Therefore, knowledge of the mass accretion rate of 
dark matter haloes is an essential ingredient in our cosmo- 
logical framework of structure formation. 

Since the collapse and virialization of dark matter 
haloes is a non-linear process, one generally resorts to nu- 
merical JV-body simulations to study the formation and 
evolution of structures in a CDM Universe. However, this 
method has a number of important drawbacks. First of all, 
numerical simulations are computationally expensive, mak- 
ing it unfeasible to explore a wide range of cosmological 
models. Secondly, computational limitations only allow sim- 
ulations of small volumes and/or a small dynamical mass 
range. This makes it virtually impossible to simultaneously 
follow the formation and evolution of objects from sub- 
galactic scale up to the scale of (super)clusters. 

An important alternative is provided by the Press- 



Schechter (PS) formalism (Press & Schechter 1974) which al- 
lows a much faster exploration of halo mass accretion histo- 
ries for a wide range of cosmologies and masses. In addition, 
the PS formalism provides a framework that allows us to 
gain insights into the physical processes involved. The main 
ansatz of the PS-formalism is to consider the initial density 
field extrapolated linearly to redshift z and smoothed on 
some typical mass scale M. Regions above a certain thresh- 
old value are then associated with collapsed objects of mass 
M at redshift z. Motivated by Birkhoff 's theorem it is there- 
fore assumed that the non-linear evolution of density pertur- 
bations, described by means of a spherical 'Top-Hat' model 
(Gunn & Gott 1972), does not influence the remainder of 
the Universe. Press & Schechter (1974) used this formalism 
to compute the mass function of dark matter haloes as func- 
tion of redshift, which has been found to be in remarkably 
good agreement with results from iV-body simulations (e.g., 
Efstathiou et al. 1988; Carlberg & Couchman 1989; Lacey 
& Cole 1994; Gelb & Bertschinger 1994; Ma 1996). 

The PS theory has been extended to give the conditional 
probabilities P(Mi,Z2\M\,z\) that a given particle at red- 
shift z\ inside a halo of mass M\ at an earlier time Z2 was 
embedded in a halo of mass M2 (Bond et al. 1991; Bower 
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1991; Lacey & Cole 1993). This extended Press-Schechter 
(hereafter EPS) formalism is easily manipulated to yield 
merger rates, halo formation/survival times, and various 
other statistics (Lacey & Cole 1993; hereafter LC93). Of par- 
ticular importance has been the construction of halo merger 
trees (e.g., Cole & Kaiser 1988; Kauffmann & White 1993; 
Somerville & Kolatt 1999; Sheth & Lemson 1999), which 
are widely used in semi-analytical models for the forma- 
tion of galaxies (e.g., Kauffmann, White & Guiderdoni 1993; 
Somerville & Primack 1999; Cole et al. 1994, 2000). 

In some cases, however, one is not interested in knowing 
the detailed distribution of halo progenitor masses, but one 
merely wants to know the rate at which the halo mass in- 
creases with time. For instance, in the disk galaxy formation 
models of Firmani & Avila-Reese (2000) and van den Bosch 
(2001a,b) the assumption is explicitly made that haloes ac- 
crete their mass smoothly; i.e., one ignores the fact that 
mass accretion involves the merging of progenitor haloes. In 
the case of disk galaxies this simplification is permitted since 
the fragility of disks suggest that mergers have not played 
an important role. In this paper we therefore use the EPS 
formalism to construct average mass accretion histories of 
dark matter haloes, which we define as the ensemble average 
(M(z)/Mo). Here M(z) is the halo mass as function of red- 
shift and Mo is the present day mass of the halo. The main 
motivation for this study is to investigate whether a simple, 
universal form exists for (M(z)/Mo). Together with the uni- 
versal profiles for the density distribution (Navarro, Frenk & 
White 1997) and angular momentum (Bullock et al. 2001) of 
CDM haloes, such universal mass accretion history provides 
a complete description of the structure and evolution of dark 
matter haloes, which can be used as a framework for detailed 
modeling of the formation of galaxies. In addition, it is to be 
expected that the accretion history of dark matter haloes is 
directly linked to the cosmic star formation history. A uni- 
versal mass accretion history might therefore proof useful in 
trying to understand the rapidly improving observations of 
the star formation rates as function of redshift. 

This paper is organized as follows. In Section ^ we start 
with a brief description of the PS-formalism and its exten- 
sion based on the excursion set formalism. Next, in Section ^ 
we describe an improved method for computing mass accre- 
tion histories, which we test against high resolution iV-body 
simulations in Section ^. In Section |s| we derive a simple, 
universal fitting formula for the average mass accretion his- 
tories of dark matter haloes, which is applicable to a wide 
range of halo masses and cosmologies. In Section ^ we dis- 
cuss a comparison of star formation and halo mass accretion 
rates, and we conclude in Section |?]. A step-by-step recipe 
for computing the average mass accretion history for a dark 
matter halo of given mass and in a given cosmology is pre- 
sented in Appendix A. 



2 THEORETICAL BACKGROUND 

In the standard model for structure formation the density 
field <5(x) = p(x)/p— 1 is considered to be a Gaussian ran- 
dom field, which is therefore completely specified by the 
power spectrum P(k). As long as 8 <C 1 the growth of the 
perturbations is linear and <5(x,t2) = <5(x, ti)D(t2)/D(ti), 
where D(t) is the linear growth factor. Once <5(x) exceeds a 



critical threshold 6° lit non-linear effects become important 
and the perturbation will start to collapse to form a virial- 
ized object (halo). In what follows we define So as the initial 
density field linearly extrapolated to the present time. In 
terms of So, regions that have collapsed to form virialized 
objects at redshift z are then associated with those regions 
for which So > 5 c {z) = S^. it /D(z)^\ 

In order to assign masses to these collapsed regions, 
the PS formalism considers the density field So smoothed 
with a spatial window function (filter) W(r;Rf), where Rf 
is a characteristic size of the filter. There is a considerable 
amount of freedom in choosing a window function, and here 
we adopt the often used spatial top-hat filter 



W(r;R f ) = 



3/(47^) 




(r < Rf) 
(r > R f ) 



(1) 



The main advantage of this filter over for example a Gaus- 
sian filter or a fc-space top-hat filter is that it is straight- 
forward to compute the mass contained within the window 
function: M — 47rp_R^/3, with p the mean mass density of 
the Universe. The ansatz of the PS formalism is that the 
probability that the density field smoothed with W(r;Rf) 
exceeds S c (z) is the same as the fraction of mass that at 
redshift z is contained in haloes with masses greater than 
M. This results in the well known unconstrained PS mass 
function for the comoving mass density of haloes: 
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(Press & Schechter 1974). Here a 2 (M) is the mass variance 
of the smoothed density field given by 



(3) 



a(M) = — I P(k) W 2 (k;R f ) k 2 dk. 



with W(k; Rf) the Fourier transform of W(r; Rf), which for 
the spatial top-hat filter used here is given by: 



W{k;R 



{kRff 



[sin(kRf) — kRf cos(kRf)] 



(4) 



The extended Press-Schechter model developed by 
Bond et al. (1991), is based on the excursion set formal- 
ism. For each point one constructs 'trajectories' S(M) of 
the linear density field at that position as function of the 
smoothing mass M. In what follows we adopt the notation 
of LC93 and use the variables S = a 2 (M) and to = 8 c (z) to 
label mass and redshift, respectively. In the limit Rf — > oo 
one has that (S, u)) = (0,0), which can be considered the 
starting point of the trajectories. Increasing 5* corresponds 
to decreasing the filter mass M, and S(S) starts to wander 
away from zero, executing a random walk. The fraction of 
matter in collapsed objects in the mass interval M, M + dM 
at redshift z is now associated with the fraction of trajec- 
tories that have their first upcrossing through the barrier 
uj = S c (z) in the interval 5*, 5* + dS, which is given by 



* Here D(z) corresponds to the linear growth factor normalized 
to unity at the present. 
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P(S,w) dS 
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(Bond et al. 1991; Bower 1991; LC93). After conversion to 
number counting, this probability function yields the PS 
mass function of equation (^) 

Since for random walks the upcrossing probabilities are 
independent of the path taken (i.e., the upcrossing is a 
Markov process) , the probability for a change AS in a time 
step Aw is simply given by equation (jH|) with S and w re- 
placed with AS and Aw, respectively. This allows one to 
immediate write down the conditional probability that a par- 
ticle in a halo of mass M 2 at z 2 was embedded in a halo of 
mass Mi at Z\ (with Zi > z<i) as 



P(Si,wi|S 2 ,w 2 ) dSi 
1 (wi 



w 2 I 



2^ (S! - S 2 )3/2 



exp 



(wi — W2) 2 



2(5i - S 2 ) 



dSi 



(6) 



Converting from mass weighting to number weighting, one 
obtains the average number of progenitors at z\ in the mass 
interval Mi , Mi + dM\ which by redshift z 2 have merged to 
form a halo of mass M 2 : 



diV 
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(M 1 ,z 1 \M 2 ,z 2 ) dMi 
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P(Si,wi|S 2 , w 2 ) 



dS 



dM 
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3 CONSTRUCTING MASS ACCRETION 
HISTORIES 

Using the EPS conditional probabilities for halo progeni- 
tor masses one can construct detailed histories of the mass 
assembly of dark matter haloes. Here we are interested in 
computing the "mass accretion histories" (hereafter MAHs) 
of dark matter haloes, which we define as *l/(Mo,z) = 
M(z)/Mo. Here M(z) is defined as the mass of the "main 
progenitor halo", and Mo is the halo mass at z = 0. 

When tracing back in time each halo breaks up in a 
number of progenitor haloes, which themselves break up in 
progenitors, etc. Given this complicated history of the mass 
evolution of a halo (often referred to as the "merger-tree"), 
we need to define what we actually mean with "the main 
progenitor" at an earlier time. We follow previous studies 
(LC93; Eisenstein & Loeb 1996; Nusser & Sheth 1999) and 
define M(z) as the main trunk of the merger tree, i.e., at 
each time step we associate M(z) with the mass of the most 
massive progenitor (hereafter MMP), and we follow that 
progenitor, and that progenitor only, further back in time. 
This way the "main progenitor halo" never actually accretes 
other haloes that are more massive than itself. Note that al- 
though at each branching point we follow the most massive 
branch, this does not necessarily imply that the main pro- 
genitor is also the most massive of all progenitors at a given 
redshift. 

In what follows we use Mo to denote the present day 
mass of a halo. Each individual time step we use M to refer 
to the parent mass for which we seek the most massive pro- 
genitor a time step Aw earlier, and we use M v to indicate 
the mass of a progenitor halo. 



3.1 Previous methods 

In an infinitesimally small time interval Aw a change AS 
results from a single merger. Therefore, an approximate 
method for constructing MAHs is to consider small time 
steps, and to assume that the change in mass associated 
with that finite time interval reflects a single merger event. 
In that case, the MMP has a mass M p > M/2 and the con- 
struction of MAHs becomes very simple: each time step one 
draws a single AS from the probability distribution 



P(AS,Aw) dAS = 



Aw 



2lF AS^ 



exp 



(Aw 2 ) 
2AS 



dAS (8) 



and one defines the mass of the main progenitor as 
max(M p , M - M p ), where a 2 (M p ) = a 2 (M) + AS. An alter- 
native approach which leads to almost similar results is to 
only accept values of AS in the range [0,cr 2 (M/2)-cr 2 (M)]. 
This method, which we call the 'binary' method, was sug- 
gested by LC93, and has been used by Eisenstein & Loeb 
(1996) to compute the minimum intrinsic scatter in the 
Tully-Fisher relation. However, as already pointed by LC93, 
this method leads to some inconsistencies regarding halo for- 
mation times, the reason for which is easily understood. The 
number density of progenitor masses diverges at small mass 
(a direct reflection of the fact that for CDM a 2 (M) —> 00 
when M — > 0), and the assumption of a single merger event 
brakes down dramatically for finite time steps, even when 
chosen very small. In other words, there is a finite, non- 
negligible, probability that the mass of the MMP is less than 
M/2. This introduces, at each time step, a systematic bias 
towards a main progenitor that is too massive, resulting in 
formation redshifts that are too high. 

An alternative method for constructing MAHs was sug- 
gested by Nusser & Sheth (1999), who draw the progenitor 
mass from the number weighted distribution function (equa- 
tion Q) with M /2 < M p < M. Although, as they show, this 
leads to better consistency with halo formation times, this 
method suffers from a similar shortcoming as again the as- 
sumption is made that the MMP is always more massive 
than M/2. However, this is only formally true in the limit 
Aw — > 0, for which the integral of dA r /dMi (equation Q) 
from M/2 to M becomes unity (i.e., it is certain that the 
MMP is more massive than M/2). 



3.2 An improved method 

The discussion above suggests that in order to construct 
MAHs using finite time steps, one has to drop the assump- 
tion of single merger events. This implies that each time step 
one needs to construct a complete set of halo masses Mi that 
are to be considered the progenitor haloes that a time step 
Aw later have merged to form the parent halo with mass M. 
An important constraint on the set Mi is that ^ . Mi = M, 
i.e, the total mass needs to be conserved at each time step. 
The MAH is then easily constructed by picking the most 
massive of Mi, and repeating the same procedure stepping 
back in time until the mass of the main progenitor is as small 
as desired. 

Unfortunately, the construction of sets of progenitors is 
not a trivial matter. In addition to conserving mass at each 
time step, a successful merger tree also has to satisfy the 
requirement that, at each time step, the distribution of the 
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Figure 1. Mass functions of progenitors (by number) for a halo with present day mass in an EdS universe with 

h = 0.65 and as = 1.0. Histograms correspond to the average number of progenitors found in 15000 random realizations of the merger 
tree. Only progenitors more massive than five percent of the parent are considered (i.e., M m ; n = 0.05M, see text), and we therefore 
only plot results for M p /Mq > 0.05. The solid lines correspond to the EPS prediction given by equation (M). Results are shown for four 
'output' times of the merger trees, the corresponding Au> and redshift of which arc indicated. In the upper panels results are shown for 
merger trees that are constructed using a fixed time step of Au> = 0.1. This results in progenitor numbers in excellent agreement with 
the direct EPS prediction. In the lower panels, progenitor masses are computed using a single time step (with Au) as indicated in the 
upper panels). In this case, too many progenitors are found compared to the EPS prediction. This indicates that for the merger tree to 
yield self-consistent results sufficiently small time steps have to be used. 



number of progenitors as function of mass is consistent with 
equation (Q). However, since the number of haloes diverges 
at very small masses, one must impose a minimum halo 
mass. Progenitors with masses below this threshold are not 
considered as individual haloes, but their mass is assumed 
to accrete smoothly onto the parent halo. As was previously 
pointed out by Somerville & Kolatt (1999; hereafter SK99), 
there seems to be no algorithm for drawing sets of progenitor 
masses that satisfies both constraints simultaneously. Kauff- 
mann & White (1993) circumvented this problem by repro- 
ducing the progenitor mass distribution exactly, but by only 
enforcing mass conservation approximately. Here we follow 
the scheme of SK99, which conserves mass exactly while 
only approximately reproducing the progenitor mass distri- 
bution. This method, termed the iV-branch method with ac- 
cretion, has been shown to yield results in good agreement 
with numerical simulations (Somerville et al. 2000). Here we 
briefly summarize the method and we refer the interested 
reader to SK99 for more details (as well as for a detailed 
discussion on other, less successful methods for constructing 
merger trees). 

The method of SK99 is based on drawing halo masses 
from the mass- weighted probability function (^) . With each 
new halo drawn it is checked whether the sum of the pro- 
genitor masses exceeds the mass of the parent M. If this 
is the case the halo is rejected and a new progenitor mass 
is drawn. Any progenitor with M v < M m i n is added to the 



mass component M acc that is considered to be accreted onto 
the parent in a smooth fashion (i.e., the formation history 
of these small mass progenitors is not followed further back 
in time). Here M m i n is a free parameter that has to be cho- 
sen sufficiently small. This procedure is repeated until the 
total mass left (M — M ac c — ^ M p ) is less than M m i n . This 
remaining mass is assigned to M acc and one moves on to the 
next time step. 

In principle, since the upcrossing of trajectories through 
a boundary is a Markov process, the statistics of progenitor 
masses should be independent of the time steps taken. How- 
ever, the SK99 algorithm is based on the single halo proba- 
bility (equation [pf), which does not contain any information 
about the set of progenitors that make up the mass of M 
(mass conservation is enforced 'by hand', by rejecting pro- 
genitor masses that overflow the mass budget). Therefore it 
is not clear whether the results are time step independent, 
and this has to be tested. In the upper panels of Figure [l] we 
plot the number distribution of progenitor masses at vari- 
ous redshifts obtained using the SK99 scheme with a fixed 
time step of /Auj = 0.1. Results are plotted for a halo with 
M = 10 12 /i _1 M in an Einstein-de Sitter (EdS) universe 
with fio = 1, £7a = 0, as = TO, and h = 0.65. Histograms 
correspond to the distributions obtained from 15000 random 
realizations of the merger tree, with M m i n = 0.05A/. Solid 
lines correspond to the EPS prediction of equation (Q). In 
the lower panels we show the same results, except that we 
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log[l+z] log[Aw] 

Figure 2. The panel on the left plots 25 random MAHs for a halo with Mo = 5.0 X 10 n /i -1 Mq in an EdS universe with h = 0.65 
and erg = 1.0. The thick solid dots with errorbars correspond to the average and standard deviations as determined from 1000 such 
realizations. A time step of Aw = 0.1 is used. The panel on the right plots curves of constant (if) (as indicated) as function of the time 
step Aw used in the construction of the MAHs. Each (ty) is determined from 1000 random realizations of the mass accretion history. For 
Aw < 0.3 the average MAHs converge to a robust result, while for larger time steps (']/) depends on the actual value of Aw used. This 
again suggests that a successful construction of MAHs requires a sufficiently small time step and we adopt Aw = 0.1 throughout. 



have now computed the progenitor masses in a single time 
step (as indicated in the upper panels); i.e., in the panels 
on the right, the upper panel shows the results when us- 
ing 100 time steps of Aw — 0.1, each time computing the 
progenitors of all previous progenitors etc. The lower panel, 
on the other hand, computed the progenitors using a single 
time step of Aw = 10.0. As can be seen, when using small 
enough time steps the number density of progenitor masses 
is in excellent agreement with equation ([?). However, when 
too large time steps are used, the method over-predicts the 
number of progenitors quite dramatically. We thus conclude, 
confirming the results of SK99, that as long as small enough 
time steps are used, the algorithm outlined above provides 
an accurate method for constructing merger trees. 

Based on this scheme we use the following algorithm to 
construct MAHs: 

(1) Choose a present day halo mass Mo and set M = Mo 
and z = 0. 

(2) Set Mieft = M and compute the progenitor redshift z p 
from Aw — 8 c (z p ) — S c (z). 

(3) Draw AS from the probability distribution (|^) and 
compute the corresponding progenitor mass M v from 
a 2 (M p ) =a 2 (M) + AS. 

(4) If M p > Mieft the progenitor mass is too big: goto 3 

(5) Set the mass of the MMP to M M mp = max[M p , M M mp] 
and the mass left in the set to Mi e ft = Mieft — M p . 

(6) If Mmmp > Mieft then we have found the MMP. In 
that case we proceed to the next time step: we set M(z p ) = 
Mmmp, M = Mmmp, z = z p , and goto 2. 

(7) Goto 3 

This procedure is repeated until the mass of the main pro- 
genitor is as small as desired. Note that, as is evident from 



step 5, we do not need to construct an entire set of progen- 
itors; we can stop once the most massive of the progenitors 
already drawn is larger than the mass left in the set, and we 
thus do not need to define a minimum progenitor mass M m i u . 
Throughout we use a fitting function for <j 2 (M) which is ac- 
curate to better than 0.5 percent over the entire mass range 
lO 6 /^ 1 M < M < 10 16 /i _1 M Q (see Appendix A). The 
power spectrum P(k) is characterized by the shape param- 
eter r. Unless stated otherwise we use the form suggested 
by Sugiyama (1995) 

r = Q h exp [-fif,(l + V2h/Qo)\ (9) 

and we adopt a baryon mass density of fif, = 0.019/i -2 
(Tytler et al. 1999). 

We now define the average mass accretion history (here- 
after AMAH) of a halo of mass Mq as 



*(Mo,*)> = -^*,(Mo,ar) 



(10) 



where the summation is over an ensemble of N random real- 
izations ^(Mo, z). Unless stated otherwise we use N = 1000 
and Aa; = 0.1. In the left panel of Figure ^ we show an exam- 
ple. The thin lines are 25 random realizations for the MAH 
of a halo with Mo = 5 x lO 11 /), -1 Mq in an EdS universe with 
as = 1.0 and h — 0.65. The solid dots with errorbars corre- 
spond to the AMAH and its standard deviation (averaged 
over 1000 random realizations). 

The right panel of Figure plots the average redshifts at 
which (*) = (0.01,0.1,0.3,0.6,0.9) as function of the time 
step Aui. Note how only for sufficiently small time steps 
(Aw Ss 0.3) the values of (^) are not time step dependent. 
This reflects the problem outlined above, that if the time 
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step is too big, the distribution of progenitor masses is no 
longer consistent with the EPS prediction of equation (Q). 
Results for other halo masses and cosmologies are similar, 
and we therefore adopt Au) = 0.1 throughout. 



4 COMPARISON WITH NUMERICAL 
SIMULATIONS 

Although the (extended) PS formalism is a valuable tool in 
the study of hierarchical structure formation, it is based on 
a number of questionable assumptions (e.g., spherical col- 
lapse). Therefore, it is essential that any statistic extracted 
using the PS formalism is tested against numerical simula- 
tions. 

4.1 Description of the simulations 

In order to test our algorithm for constructing MAHs we 
compare our results with the high resolution ACDM simu- 
lation of the GIF project (Kauffmann et al. 1999; Diaferio 
et al. 2001). This simulation is performed using the par- 
allel adaptive particle-particle particle-mesh (AP 3 M) code 
called Hydra (Couchman, Thomas & Pearce 1995) and fol- 
lows the evolution of 256 3 (~ 16.8 x 10 6 ) particles in a 
141.3/i _1 Mpc box. The cosmological parameters used are 
Slo = 0.3, Qjv = 0.7, h = 0.7, and F = 0.21. The power- 
spectrum of density fluctuations is normalized to a$ = 0.9, 
in agreement with the observed abundance of massive clus- 
ters (Eke, Cole & Frenk 1996). Particle mass and softening 
length are 1.4 x 10 10 ft _1 M and 30/i _1 kpc, respectively. 

Particle positions and velocities are stored at 43 differ- 
ent redshifts between z = 12 and z = 0. At each output time, 
haloes are identified using the standard friends-of-friends 
(FOF) percolation algorithm (Davis et al. 1985), where a 
linking length of 0.2 times the mean inter-particle separation 
is adopted. The mass associated with each of these haloes, 
denoted by Mfof, is simply given by the number of particles 
linked together times the single particle mass. 

A halo at redshift 22 is defined as a progenitor of a halo 
at redshift z\ < 22 if (i) more than half its particles are 
included in the halo at zi, and (ii) its most bound particle 
is also included. Using lists of progenitors for all haloes and 
for all output timesjf] a MAH is constructed as follows. We 
start with a halo at z=0 and identify all its progenitors at 
the previous time step z\. The new halo mass M(z\) is then 
simply defined as the mass of the most massive of these pro- 
genitors, and this procedure is repeated for all subsequent 
output times. A more detailed description of the construc- 
tion of the halo merger trees from these GIF simulations can 
be found in Kauffmann et al. (1999). 




log[l+z] 

Figure 3. The thin lines represent mass accretion histories for 
haloes in the ACDM GIF simulation. In upper panel plots the 
MAHs for 25 (randomly selected) haloes with 2.0 X 10 14 /i _1 M Q < 
M FO f < 4.0 X 10 14 h. _1 Mq. The thick lines and dots with er- 
rorbars correspond to the average MAH computed using our 
method based on the EPS formalism for a halo mass of Mq = 
2.75 X 10 14 /i -1 Mq, corresponding to the average FOF mass of 
the 25 haloes in the simulation. As can be seen, haloes in the 
simulation seem to form a little bit earlier as predicted by the 
EPS formalism. In the lower panel we plot the MAHs of those 
FOF groups in the above mass interval which have dM/dz < 
at all times; only 6 haloes obey these criteria. As can be seen, the 
MAHs of these haloes are not statistically different from those in 
the central panel, and equally inconsistent with the EPS predic- 
tion (see discussion in text). 



4.2 Comparison of mass accretion histories 

Here we compare the MAHs of dark matter haloes in the 
simulation to the AMAHs com puted using the EPS formal- 



ism as outlined in Section 3.2. 

In the upper panel of Figure [] we plot the MAHs for 
25 (randomly selected) haloes in the simulation which at 



t Available at http://www. mpa-garching.mpg.de/GIF/ 



2 = have masses in the range 2.0 x W 14 h 1 Mq < Mfof < 
4.OxlO 14 ft _1 M . This mass range is chosen to ensure haloes 
with large numbers of particles, so that we are less sensitive 
to resolution issues. We only accept haloes for which the 
MAHs can be traced back to the point where < 0.01 (i.e., 
in some cases no progenitors can be identified in the sim- 
ulation before the mass of the main progenitor has fallen 
below one percent of the present day mass). The thick solid 
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dots with errorbars correspond to the AMAH computed 
using the EPS formalism (for exactly the same cosmolog- 
ical parameters as used in the simulation) for a halo with 
Mo — 2.75 x 10 14 /i _1 Mq which corresponds to the aver- 
age mass of the 25 haloes in the simulation. Although the 
overall shape of the MAHs and the amounts of scatter are 
fairly similar, the MAHs of the haloes in the numerical sim- 
ulation are systematically offset to higher redshifts (i.e., the 
haloes in the simulation seem, on average, to form somewhat 
earlier) . 

One possible explanation for this inconsistency is that 
the haloes in the simulation often have progenitors that are 
more massive than the parent. This can have several causes. 
Haloes in the simulations are susceptible to tidal stripping, 
which can cause the mass of the main progenitor to actu- 
ally decrease with time. In addition, during a (high velocity) 
encounter two haloes may temporarily be linked together 
by the FOF algorithm, which causes M(z) to increase and 
subsequently decrease again. These effects are not modeled 
by the EPS formalism, which only allows halo masses to 
grow with time. We can investigate the effect this has on 
the statistics of our MAHs by only selecting haloes from the 
simulation that, at each time step, have dM/dz < 0. Only 
6 of the 35 haloes in our mass interval obey this criterion, 
and their MAHs are plotted in the lower panel of Figure pi 
Although the number statistics are poor, it is evident that 
these MAHs are not statistically different from those that 
occasionally show dM/dz > 0, and they show a similar in- 
consistency with respect to the expected AMAH. We thus 
conclude that the inconsistencies between the MAHs from 
the simulations and those from the EPS formalism are not 
related to the fact that the latter does not take possible mass 
loss into account. A more likely cause for the inconsistency 
is discussed in Section 1.4. 



4.3 Comparison of halo formation times 

Another useful statistic for comparison is the distribution of 
halo formation redshifts. We use the definition of LC93 and 
define the halo formation redshift, Zf, as the redshift where 
$ = 0.5, i.e., where the mass of the main progenitor is half 
the present day mass. 

Lacey & Cole (1993) presented two methods for com- 
puting the distribution of formation times. The first is based 
on the number weighted mass distribution of progenitor 
haloes (equation J?)). As argued by LC93, integrating this 
equation from Mo/2 to Ado gives the probability that the 
progenitor mass is more massive than Mo/2, which is equal 
to the probability that the halo formation time was earlier 
than this. Upon defining the scaled variables 



S = 



and 



q 2 (M) -a 2 (Mo) 
r 2 (M /2) -a 2 (M ) 



S c (zf)-S c (0) 
y/a 2 (Mo/2)-a 2 (Mo) 



(11) 



(12) 



one can write the probability distribution for halo formation 
times as: 



where M is solved from equation (]ll|). The advantage of 
using the variables S and w/ is that P(u)f) depends only very 
mildly on mass and cosmology (this dependence is largely 
absorbed by the variables themselves). 

Another method for computing halo formation times 
is to use the actual MAHs themselves, and to identify the 
redshift at which the main progenitor mass equals half the 
present day mass. LC93 used their 'binary' Monte- Carl o 
method for constructing MAHs (discussed in Section 3.1), 



/2tt 



Mo 
M 



SB/2 g3/2 







j exp 




25. 



dS, (13) 



and found a distribution P(u>f) that was offset from that of 
equation ( |l3| ) to higher values of uif. Lacey & Cole (1994) 
determined the distributions of halo formation times of dark 
matter haloes in simulations with scale-free power-spectra, 
and found them to be in excellent agreement with the an- 
alytical prediction of equation (|l^), but inconsistent with 
P(Cbf) derived using the Monte-Carlo method. As indicated 
in Section [Tl], the 'binary' Monte-Carlo method used by 
LC93 is based on the false premise that the most massive 
progenitor has mass M v > M/2, which leads to a systematic 
offset of halo formation redshifts to too high values. 

We now re-examine the issue of halo formation times us- 
ing the high-resolution ACDM simulation and our improved 
method for constructing MAHs. In order to improve accu- 
racy we use linear interpolation between the time steps that 
bracket $ = 0.5 (for both the haloes in the simulations as 
well as for the MAHs constructed using the EPS formalism). 
We convert zj to the scaled variable uif, so that P(u>f) de- 
pends only very weakly on cosmology and halo mass. This 
allows us to compare distributions of uif for relatively large 
mass intervals, providing better statistics. 

We compute the distribution of uif for two samples of 
haloes. The first consists of 5044 haloes with present day 
masses 5.6 x 10 n /i _1 M < M PO p < 1.12 x W 12 ^ 1 M ). 
This mass range corresponds to haloes that, at z = 0, con- 
sist of between 40 and 80 particles. The resulting distri- 
bution of uif is indicated by the hatched histogram in the 
left panel of Figure ^. The solid curve corresponds to the 
distribution of uif derived from 2 x 10 5 MAHs for a halo 
with Mo = 8.4 x 10 h~ M (corresponding to the average 
mass of the 5044 haloes in the simulation), using the same 
cosmological parameters as in the simulation. The two dis- 
tributions are in good agreement, although the mean for the 
simulated haloes is slightly offset to higher uif. The dashed 
curve corresponds to P(uif) of equation (|l3|), also for a mass 
of 8.4 x lO 11 /!" 1 M . This distribution is in excellent agree- 
ment with that derived using the MAHs. This shows that our 
method for constructing MAHs is, within the EPS frame- 
work, self-consistent, and does not suffer from the inconsis- 
tencies found when using the LC93 'binary' method. 

The right panel of Figure ^ plots similar results but now 
for the 6042 haloes in the mass range 2.0 x lO 11 ^ 1 M < 
Mfof < 4.0 x 1O 12 /i _1 M . The agreement with the MAHs is 
somewhat poorer as in the case of the lower mass range, with 
significantly higher formation redshifts for the haloes in the 
simulation. This is, off course, another manifestation of the 
inconsistencies found between the MAHs in the simulations 
and those computed using the EPS formalism (cf., Figure ^). 

4.4 Ellipsoidal collapse 

In addition to the inconsistencies with the MAHs indi- 
cated above, various studies in the past have pointed out 
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Figure 4. The hatched histograms correspond to the distribution of halo formation times (parameterized by the scaled variable u>t) for 
haloes in the simulation (with FOF masses as indicated in the two panels). The solid lines are the distributions of uit computed from 
2 X 10 5 random realizations of the MAH for a halo with M = 8.4 X lO^h' 1 M Q (left panel) and M = 1.1 X 10 13 h _1 M (right panel}. 
These masses correspond to the average halo masses in the two simulation samples. Dashed lines correspond to P(£jf) of equation ( |13|) 
for the same masses, and are in excellent agreement with the distributions obtained directly from the MAHs. However, the formation 
times in the simulation are offset from those of the EPS formalism, whereby the discrepancy is larger for more massive haloes. 



that the standard, unconditional PS mass function (equa- 
tion El) over (under) predicts the number of low (high) 
mass haloes when compared to numerical simulations (e.g., 
Jain & Bertschinger 1994; Tormen 1998; Gross et al. 1998; 
Governato et al. 1999). This is generally interpreted as a 
consequence of the assumption of spherical collapse, and 
numerous studies have shown that considering ellipsoidal 
rather than spherical collapse brings the PS mass function 
in much better agreement with simulations (e.g., Monaco 
1995; Bond & Myers 1996; Audit, Teyssier & Alimi 1997; 
Lee & Shandarin 1998; Sheth & Tormen 1999; Lanzoni, Ma- 
mon & Guiderdoni 2000; Sheth, Mo & Tormen 2001; Jenkins 
et al. 2001). 

We can investigate a modification of the spherical col- 
lapse model, by multiplying the critical collapse density S c 
(as defined in Section H) with a fudge factor a and by exam- 
ining how a depends on halo mass, if at all. One of the ad- 
vantages of using the scaled variable £>/ is that P(u>f) is inde- 
pendent of a, at least for the MAHs computed using the EPS 
formalism. In the case of the haloes in the simulation one has 
that u>f oc a, and one can therefore immediately re-scale the 
histograms in Figure ^ to determine the best-fit value of a. 
Doing so we find a ~ 0.94 and a ~ 0.82 for the mass intervals 
plotted in the panels on the left (M = 8.4 X lO^h' 1 M ) 
and right (Mo = 1.1 x 10 13 /i -1 M Q ), respectively. Similarly, 
we find that for a ~ 0.8 the AMAH plotted in Figure | 
(M = 2.75 x 10 14 /i _1 M Q ) is in excellent agreement with 
the MAHs found in the simulation. Thus it seems that a 
modification of the spherical collapse model whereby S c de- 
creases with increasing halo mass can bring the EPS MAHs 
in excellent agreement with the simulations. Interestingly, 
as shown by Sheth, Mo & Tormen (2001), this is exactly the 
kind of behavior one expects if one takes into account that 



haloes are ellipsoidal rather than spherical. These authors 
obtain a modified critical collapse overdensity given by 



5 ce {M,z) = 5 c (z) 1 + 0.47 



\M) 



(14) 



Here 8 c (z) is the standard value for the spherical collapse 
model. This modification results in halo mass functions that 
are in excellent agreement with those found in simulations 
(Sheth & Tormen 1999; Jenkins et al. 2001). 

Thus, two separate statistics, the unconditional halo 
mass function and the halo formation times (i.e., the MAHs), 
both suggest a critical collapse density that increases with 
decreasing halo mass. Therefore it might be worthwhile to 
try and incorporate a mass-dependent critical collapse den- 
sity in the EPS formalism. This requires, what pundits call, 
determining the upcrossing statistics for a moving barrier. 
However, as discussed by Sheth & Tormen (2001) no analyt- 
ical expression for the conditional probability function (i.e., 
a moving barrier equivalent of equation |J) is known for a 
barrier of the form of equation (|l4|), and one has to resort 
to either an approximate fitting function, or one has to use 
time-consuming Monte-Carlo simulations to determine the 
upcrossing statistics. Unfortunately, as discussed in detail 
by Sheth & Tormen (2001), neither of these two methods 
are appropriate for constructing merger trees. Therefore, in 
what follows we adhere to the standard, spherical collapse 
model, but we caution the reader that, taking the numerical 
simulations at face value, the MAHs thus derived contain 
slight, mass-dependent inaccuracies. 
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log[l+z] 

Figure 5. Average MAHs for various cosmologies (all with erg = 1.0 and h = 0.65). Results are shown for five masses: 5.0 X I0 9 h~ 1 M 
(open triangles), 5.0 X 10 10 /t _1 M Q (open squares), 5.0 X lO 11 /^ 1 Mq (crosses), 5.0 X 10 12 /i _1 M (tripods), 5.0 X 10 13 /i _1 M (open 
circles). Averages are determined from 10 3 random MAH realizations. Solid lines correspond to the universal MAH computed using the 
recipe given in Appendix A, and provide an excellent fit to the average MAHs. 



5 A UNIVERSAL FORM FOR THE MASS 
ACCRETION HISTORIES 

Figure |sj plots the average MAHs for various halo masses 
(different symbols) and cosmologies (different panels). Up- 
per panels correspond to ACDM cosmologies with fio+fiA = 
1, while lower panels are for OCDM cosmologies without 
cosmological constant (in both cases we adopt h = 0.65 and 
erg = 1.0). These plots clearly show the well-known behav- 
ior that smaller mass haloes form earlier, which is a direct 
reflection of the fact that a(M) increases with decreasing 
mass. The cosmology- dependence of the MAHs is easily un- 
derstood if one takes into account how the mass variance 
<t(M) and the linear growth factor D(z) depend on cosmol- 
ogy. On the mass scales of interest, a decrease in fio causes 
a(M) to decrease (for fixed as), which implies lower for- 
mation redshifts (i.e., slower accretion rates). At the same 
time, a decrease in fio causes an increase in D(z) (at fixed 
redshift), so that a time step Aui implies a larger Az. This 
drives the MAHs to higher formation times. The net result of 
a decrease in fio, therefore, depends on which of these two ef- 
fects dominates. Since Aa/AM decreases with fio the a(M)- 
effect on z/ is stronger for less massive systems. Therefore, 
one expects the a(M)-effect to dominate for small enough 
masses, resulting in a decrease of Zf. Furthermore, the in- 
crease of D(z) with decreasing fio is stronger for a Universe 
with £Ia = than for one with Qa = 1 — fio, so that the 
mass scale below which Zf decreases with decreasing fio is 
higher in an open cosmology compared to a A-cosmology. 
This behavior is nicely reproduced by the MAHs plotted in 
Figure |^. In the upper panels, decreasing fio only very mildly 
affects the MAH of a 5 x 10 13 /i _1 M n halo. For this mass 



the two effects mentioned above largely cancel each other. 
For less massive haloes the <r(M)-effect dominates, causing 
the haloes to form later in lower-fio cosmologies. In fact, for 
fio =0.1 the mass dependence of MAHs is much reduced 
compared to the EdS cosmology. In OCDM cosmologies, the 
D(z)-effect is relatively stronger and now it are the low mass 
systems that are hardly affected, while more massive haloes 
increase their formation redshift with decreasing fio. 

After experimenting with a variety of fitting functions, 
we find that the AMAHs are well fitted by the following 
simple form: 



log<*(M ,z)) = -0.301 



log(l + 3) 



log(l + z/) 



(15) 



where Zf and v are free fitting parameters. Note that with 
this definition Zf corresponds to the formation redshift as 
defined in Section 4.S (i.e., M(zf) = Mo/2). In what follows 
we shall refer to equation ( JTHj ) as the 'universal' mass accre- 
tion history. Recently Wechsler et al. (2001) investigated the 
mass accretion histories of individual dark matter haloes in 
a ACDM simulation which they fitted with *(M , z) = e~ az 
(with a a free fitting parameter). In Appendix B we compare 
this one-parameter fitting function to the universal MAH of 
equation (^), and we show how a may be estimated from 
the formation redshift Zf. 

In order to investigate how the scale-parameters v and 
Zf depend on mass and cosmology we proceed as follows: we 
randomly draw values for fio, Mo, and as from the inter- 
vals 0.1 < fio < 1.0, 10 9 /i _1 M < M < 10 14 ft _1 M and 
0.5 < as < 1.5. For each of these models we compute the 
AMAH from 10 3 random realizations of the MAH, and we 
find the best-fit values of v and Zf. In total we construct 50 
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Figure 6. The correlation between the best-fit values of the scale- 
parameters v and zj of the universal MAH, obtained by fitting 
equation Jl5[ ) to the AMAHs of dark matter haloes with a large 
range of masses and for a wide variety of cosmologies. Errorbars 
correspond to the formal errors returned by the fitting routine 
used. Note that v and log(l + Zf) are strongly correlated suggest- 
ing that a single parameter suffices to describe the average mass 
accretion history of a dark matter halo (see Appendix A). 



AMAHs with fl A = and another 50 with fi A = 1 - fio- As 
shown in Figure |(| the best-fit values of v and Zf are strongly 
correlated, suggesting that a single parameter suffices to de- 
scribe the average MAHs. In Appendix A we present a step- 
by-step recipe for computing v and 2/ as function of halo 
mass and cosmology directly. 

The solid lines in Figure [| correspond to the univer- 
sal mass accretion histories with v and 2/ computed using 
this recipe. Only in the extreme case with Qo = 0.1 and 
JIa =0.0 (lower right panel), does the universal MAH fail 
to accurately fit the average MAHs. In all other cases, the 
universal MAH is in excellent agreement with ave rage MAHs 
obtained using the method described in Section |3.2|. 



6 APPLICATION: THE ACCRETION RATE 
OF DARK MATTER HALOES 

Using the universal mass accretion history, the average rate 
at which haloes of mass Mo accrete mass can be written as 



AM 



, , d* 
~dt~ ' 



(16) 



We can use this to compute the total comoving dark mat- 
ter accretion rate of all haloes that at z = have masses 
between Mi and M%, by weighting each halo by the present 
day, comoving number density n(M, z — 0): 



^(*;Mi,M a ) = 



M 2 



Mi 



~A~t ^ 



dn 
dlnM 



(M,z = 0) dM (17) 



Here p(z) is defined as the comoving mass density at red- 
shift z of all main progenitor haloes that by z = have 
evolved to become haloes with Mi < Mo < Mi. Using the 
universal MAH of equation (Jl5|) and the PS mass function 
of equation (g) , this integral is easily computed numerically. 

If we make the simplifying assumption that all baryons 
inside haloes with present day masses in the interval Mi < 
Mo < M2 were instantaneously turned into stars the mo- 
ment they were accreted, multiplying equation ( |l7| ) with 
the universal baryon fraction /b ar gives the comoving star 
formation rate dp*/dt as function of redshift. For Mo <; 
10 13 /i _1 Mq the cooling time is longer than the Hubble time, 
and such systems are therefore not expected to contribute 
significantly to the cosmic star formation rate. Haloes with 
M <, 10 10 /i _1 M have virial velocities Kir & 30 kms" 1 , 
and a typical background UV radiation field can prevent the 
gas from cooling (e.g., Babul & Rees 1992; Kepner, Babul 
& Spergel 1997). Therefore it is to be expected that the ma- 
jority of star formation occurs in haloes in this mass range. 
We therefore set Mi = lO 10 /^ 1 M and M 2 = 10 13 /t _1 M , 
and compute dp*/di using /b ar = 0.019^ h~ 2 (Tytler et 
al. 1999). 

The results are shown in Figure (7] for three different cos- 
mologies, all with h — 0.65 and as = 1.0 (solid lines). The 
short-dashed, dotted, and long-dashed curves plot the sepa- 
rate contributions from the mass ranges 10 < log(Afo) < 11, 
11 < log(Mo) < 12, and 12 < log(M ) < 13, respectively. 
In all three cosmologies dp*/di increases rapidly at low red- 
shift, peaks in the redshift interval 1 z ^ 3, and then de- 
clines (the steepness of which depends on cosmology). This 
is in good agreement with observations of the cosmic star 
formation rate, which increases by over an order of mag- 
nitude from z = to z = 1, and which seems to peak at 
1 <, z <; 2 (e.g., Lilly et al. 1996; Madau et al. 1996; Steidel 
et al. 1996; and references therein). Clearly our assumption 
of instantaneous star formation with 100 percent efficiency 
is a severe over-simplification. In reality, there will be a de- 
lay between accretion and star formation set by the cooling 
and free fall time scales of the halo. In addition, not all 
baryons partake in star formation, as present day gas mass 
fractions in galaxies are clearly not zero. Furthermore, var- 
ious processes can temporarily quench or enhance star for- 
mation compared to the cooling rate, and feedback processes 
can cause baryons to cycle through multiple star formation 
episodes. A more elaborate comparison with the observed 
cosmic star formation rate will have to take all these effects 
into account. Nevertheless, it is reassuring that our over- 
simplified assumptions already yield results that reproduce 
the main characteristics observed. This suggests that indeed 
the baryonic mass accretion rate is the main driving force 
for the cosmic star formation history, and that the universal 
MAH derived here may proof a useful tool in modeling the 
history of star formation in the Universe. 



7 CONCLUSIONS 

We have presented an improved method for determining 
the mass accretion histories (MAHs) of dark matter haloes, 
based on the A^-branch merger-tree construction algorithm 
of Somerville & Kolatt (1999). As we have shown, this yields 
MAHs with formation times that are in excellent agreement 
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redshift 

Figure 7. The baryonic mass accretion rates defined by equation ( |l7| ) multiplied with the Universal baryon fraction /bar- Solid lines 
correspond to the comoving accretion rates integrated over all haloes with present day masses 10 10 /i _1 Mq < Mq < 10 13 /i -1 Mq. In 
addition we plot the contribution to these curves from three separate mass intervals: lO 10 /! -1 M < M < lO 11 /! -1 M (short-dashed 
lines), lO 11 /! -1 Mq < M < 10 12 /i _1 Mq (dotted lines), and lO 12 /^ 1 Mq < Mo < lO 13 /^ 1 Mq (long-dashed lines). Results are shown 
for three different cosmologies as indicated (all with h = 0.65 and ag = 1.0). The rapid increase at low z to a maximum accretion rate 
at z ~ 2 is in good agreement with observations of the cosmic star formation history (see discussion in text) . 



with direct estimates based on the extended PS formalism. 
This solves an inconsistency that hampered previous meth- 
ods for constructing MAHs, which were based on a binary 
method where the assumption was made that the most mas- 
sive progenitor is always more massive than half the mass of 
the parent. This, however, is a poor assumption and results 
in halo formation redshifts that are too high. 

The MAHs and halo formation times obtained using 
our improved method are in reasonable agreement with high 
resolution numerical simulations. The small discrepancies 
found seem to be larger for more massive haloes. Such mass 
dependence is also found when comparing the PS mass func- 
tion with simulations. Various authors have suggested that 
this mass function discrepancy can be solved by adopting an 
ellipsoidal collapse model, rather than the spherical collapse 
model used in the standard PS formalism. Interestingly, un- 
der ellipsoidal collapse one predicts the critical collapse den- 
sity to be higher for less massive systems, which is consistent 
with our mass dependence of the halo formation time dis- 
crepancies found. 

We have shown that the average MAHs of dark mat- 
ter haloes have a universal functional form. The depen- 
dence of mass and cosmology is absorbed by two parameters 
that can be computed accurately using simple fitting func- 
tions. Together with the universal density profile of CDM 
haloes (Navarro, Frenk & White 1997), and the universal 
angular momentum distribution of CDM haloes (Bullock et 
al. 2001), this universal mass accretion history provides a 
complete set of simple equations that can be used to model 
the structure and formation of the population of dark mat- 
ter haloes. The universal mass accretion history is especially 
useful in modeling the formation of disk galaxies, where 
the detailed effects of merging are not important (e.g., the 
models of Firmani & Avila-Reese 2000 and van den Bosch 
2001a,b). In addition, the universal MAH allows a straight- 
forward computation of the mass accretion rate of dark mat- 
ter haloes, which is expected to drive the cosmic star for- 



mation rate, and allows a fast but accurate investigation of 
mass and/or cosmology dependencies without the need for 
constructing ensembles of actual mass accretion histories. 
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APPENDIX A: A RECIPE FOR COMPUTING 
THE AVERAGE MASS ACCRETION HISTORY 
OF DARK MATTER HALOES 

As shown in Section |E] the average MAH of dark matter 
haloes is well fitted by the universal form: 



log(*(M ,2)) = -0.301 



log(l + z) 



Iog(l + 2/) 



(Al) 



with Zf and v two scale-parameters that depend on halo 
mass and cosmology, and which are strongly correlated (cf. 
Figure m . 



3 - 



0) 
T3 

O P 

£ 2 



1 - 







_ 1 1 1 1 1 1 


i ■ ■ ■ ■ i ■ ■ ■ ■ _ 

: 




o 

a _ 


— 


So 

y - 


1 1 1 1 1 1 


, , , i , , , , i , , , , ■ 







2 

z f (fit) 



Figure Al. The values of Zf derived from equation (A6) with 
/ = 0.254 ver sus the best-fit values of zj determined from fitting 
equation (Al) to the average MAHs. Results are shown for 100 
haloes that differ in both mass and cosmology (see Section pi). 



In order to compute Zf directly from EPS theory we 
consider the cumulative probability distribution of single 
trajectory formation redshifts given by LC93: 



P(z > 2/ 1 Mo) = erfc 



5 c (z f ) - 6c(Q) 



v / 2{a 2 {fM )-a 2 (M ) 



(A2) 



with erfc(:r) the complementary error function and / = 0.5. 
Here 5 c (z) — S^ lit [Q,(z)]/D(z) with D(z) the linear growth 
factor normalized to unity at z — (see Peebles 1980) and 
<5° rit is the critical threshold for spherical collapse, which is 
well approximated by 



S° lit [Q(z)] 
with 



0.15(127r) 2/3 n(z) 



0.0185 if fi < 1 and Qa = 
0.0055 if Q + Oa = 1 



(e.g., Navarro, Frenk & White 1997), and 

n (i + z) 3 



fi(z) = 



fi A + (1 - fio - + z) 2 + fi (l + z) 



(A3) 



(A4) 



(A5) 



The median value for Zf of equation (A2) is easily obtained 
by solving for the root of 



8 c (z f ) = Sc(0) + 0A77 v / 2[a 2 (fM o ) ~ a 2 (M Q )] 



(A6) 



As pointed out by LC93, this is not the same as the median 
halo formation redshift, as it does not nece ssar ily follow the 
main progenitor. However, since equation (A6) at least con- 



tains the proper scaling with cosmological parameters, one 
might hope to be able to use this simple equation to model 
the best- fit values of Zf by tuning the parameter /. We find 
excellent agreement with the best-fit values of 2/ (deter- 
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Figure A2. The left panel plots the ratio of Xexp an d Xuni as function of the best- fit value of zj for 100 haloes with different masses 
and in different cosm olog ies (see Section |jj). As can be seen, on average Xcxp > Xuni f° r z f < 1-5, which implies that the Universal fitting 
function o f eq uation o provides a better fit to the AMAH of CDM haloes that form relatively late, whereas the exponential MAH of 
equation fell) typically provides a better fit in cases with Zf > 1.5. In the panel on the right we plot the best-fi t va lues of a versus the 
best-fit values of zt, which are strongly correlated. The solid line corresponds to the fitting function of equation (B2j), which can be used 
to convert Zf computed using the recipe in Appendix A to a. 



mined by fitting equation (Al) to the average MAHs) for 
/ = 0.254 (see Figure O). 



Motivated by the strong correlation between the best- 
fit values of v and Zf we find that the power-index v is well 
approximated by 

v = 1.211 + 1.858 log[l + z f ] + 0.308 

0.032 logbMoAlO 11 /! -1 M )] (A7) 

(rms error of 9.0 x 10" 4 ). Equations (|X|) (with / = 0.254) 
and (A7) give the scale-parameters Zf and v for any halo 



mass and cosmology and therewith completely specify the 
universal mass accretion history (Al). 

Whenever we compute a(M) we use the following fitting 
function: 



ct(M) = erg 



m 

/(us)' 



(A8) 



Here us = 32r (with T the power spectrum shape parame- 
ter), 



4 / M\ 1/3 
„ = 3.804 x 10- r(-) , 

with M in units of h~ x M©, and 
/(it) = 64.087 [l + 1.074 u ' 3 - 1.581 u A - 
0.954 u °' 5 - 0.185 u 06 ]- 10 



(A9) 



(A10) 



This fitting function, which is accurate to better than 
0.5 percent over the mass range 10 6 /i _1 M Q < M < 
10 16 /i _1 Mq, is only valid for a spatial top hat filter, and 
based on a power spectrum P(k) = kT 2 (k) (i.e., we assume 
that the initial power spectrum has a Harrison-Zeldovich 
form P(k) oc k) with T(k) the transfer function given by 



Bardeen et al. (1986): 

Ml + 2.34,) 
v ' 2.34(7 

[l + 3.89g + (16.1g) 2 + (5.46q) 3 + (6.71g) 4 ] " 1/4 (All) 
Here q = k/Y, with k in ftMpc -1 . 



APPENDIX B: AN ALTERNATIVE FORM FOR 
THE UNIVERSAL MASS ACCRETION 
HISTORY 

After this paper was submitted Wechsler et al. (2001; here- 
after WBPKD) presented a similar investigation into the 
mass accretion histories of CDM haloes. They used the ex- 
ponential form 



V(M ,z) = e 



(Bl) 



to fit the mass accretion histories of individual haloes in 
a ACDM simulation. There are subtle differences between 
the fitting function suggested here and the one used by 
WBPKD, warranting a closer inspection which, if any, pro- 
vides a better fit. To that extent we have fitted the AMAHs 
of the 100 haloes (with different masses and different cos- 
mologies) presented in Section |E| with the exponential MAH 
of equation (Bl). In the left panel of Figure A2 we plot 



log(Xoxp/Xunijas function of z f . Here Xox P and Xuni cor- 
respond to the best-fit values of % 2 f° r t ne fitting func- 
tions of equations (Bl) and (^), respectively. As can be 
seen, log(Xex P /Xuni) an d Zf are fairly strongly correlated to 
the extent that haloes that form relatively late (zf 1.5) 
are on average better fit by the Universal MAH of equa- 
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tion (|l5|), whereas the opposite is true for haloes that form 
early (z/ J; 1.5). 



In the right panel of Figure AS we plot the best-fit val- 
ues of a as function of the best-fit values of Zf. As can be 
seen there is a strong correlation between these two param- 
eters, which is well fitted by 



1 1.43 



(B2) 



(solid line). The fact that a and zt are so well correlated 



means that one can use equation (B2) and the recipe in Ap- 
pendix A to estimate the value of a for the average MAH of 
a dark matter halo of arbitrary mass and cosmology. This 
can, for example, be used to improve accuracy by using equa- 



tion (Bl) instead of ( |15| ) as an analytical description for the 
AMAH of haloes with Zf <; 1.5. 

Another advantage of having a direct way to compute 
a from the recipe outlined in Appendix A is that WBPKD 
have shown that a is strongly correlated with the concen- 
tration of the final halo; more concentrated haloes form on 
average earlier. WBPKD have shown that a good fit is ob- 
tained with Cvir = 8.2/ a, where c v ir = r v i r /r a with r v i r the 
virial radius of the halo, and r s the characteristic radius of 
the NFW (Navarro, Frenk & White 1997) halo density pro- 



file. Using the recipe in Appendix A and equation (B2) one 
can thus compute the average concentration of a dark mat- 
ter halo of any mass and for any cosmology. There are two 
small caveats here. First of all, the relation between c v j r and 
a has only been tested for one cosmology (with f2o = 0.3, 
Qa = 0.7, h — 0.7 and erg = 1.0) and it remains to be 
seen whether this also holds for other cosmologies. However, 
since the results of WBPKD strongly suggest that the halo 
concentration is set entirely by its MAH (i.e., the scatter in 
Cvir for haloes of a fixed mass is consistent with the scatter 
in MAHs), it seems likely that this will be the case. Sec- 
ondly, the relation between zj and a of equation ( |B2| ) is 
based on EPS MAHs, whereas the relation between c v ir and 
a is based on MAHs fitted directly to simulations. Since 
the EPS MAHs and those extracted from simulations show 
some inconsistencies (see Section ^) there will be a system- 
atic, but small, error in the values of c v ir thus derived. Since 
the discrepancy between EPS and simulations depends on 
mass and cosmology, the same applies for the amplitudes 
of these errors. Nevertheless, since for most cases the error 
will be relatively small compared to the scatter in c v ir, this 
method allows one to compute halo concentrations for ar- 
bitrary halo mass and cosmology to sufficient accuracy for 
most purposes. 



